#include "pbc.h"
#include "const.h"
#include "xyz.h"

namespace kuai {

	PBC::PBC() {
		ax = by = cz = LARGE_NUMBER;
		bx = bz = cx = 0;
	}

	PBC::PBC(const RealNumber& a, const RealNumber& b, const RealNumber& c) 
	{
		ax = a;
		by = b;
		cz = c;
		bx = bz = cx = 0;
	}

	PBC::PBC(const RealNumber& a, const RealNumber& b, const RealNumber& c, const RealNumber& alpha, const RealNumber& beta, const RealNumber& gamma)
	{
		RealNumber cosA = cos(deg2rad(alpha)), cosB = cos(deg2rad(beta)), cosG = cos(deg2rad(gamma));
		RealNumber sinB = sin(deg2rad(beta));
		ax = a;
		bx = b * cosG;
		cx = c * cosB;
		by = b * sqrt(1-cosA*cosA -cosB*cosB - cosG*cosG + 2*cosA*cosB*cosG) / sinB;
		bz = b * (cosA - cosB*cosG) / sinB;
		cz = c * sinB;
	}

	RealNumber PBC::volumn() const {
		return ax * by * cz;
	}
	
	XYZ PBC::image(const XYZ& v) const {
		XYZ result(v);
		update(result);
		return result;
	}
	void PBC::update(XYZ& v) const {
		while (v.y < -by/2) {
			v.x += bx;
			v.y += by;
			v.z += bz;
		}
		while (v.y > by/2) {
			v.x -= bx;
			v.y -= by;
			v.z -= bz;
		}
		while (v.z < -cz/2) {
			v.x += cx;
			v.z += cz;
		}
		while (v.z > cz/2) {
			v.x -= cx;
			v.z -= cz;
		}
		while (v.x < -ax/2) {
			v.x += ax;
		}
		while (v.x > ax/2) {
			v.x -= ax;
		}		
	}
	void PBC::update(RealNumber v[]) const {
		update(reinterpret_cast<XYZ&>(v));
	}
	
	void PBC::shift(const XYZ& v, int result[]) const {
		XYZ vv(v);
		result[0] = result[1] = result[2] = 0;   
		while (vv.y < -by/2) {
			result[1] += 1;
			vv.x += bx;
			vv.y += by;
			vv.z += bz;
		}
		while (vv.y > by/2) {
			result[1] -= 1;
			vv.x -= bx;
			vv.y -= by;
			vv.z -= bz;
		}
		while (vv.z < -cz/2) {
			result[2] += 1;
			vv.x += cx;
			vv.z += cz;
		}
		while (vv.z > cz/2) {
			result[2] -= 1;
			vv.x -= cx;
			vv.z -= cz;
		}
		while (vv.x < -ax/2) {
			result[0] += 1;
			vv.x += ax;
		}
		while (vv.x > ax/2) {
			result[0] -= 1;
			vv.x -= ax;
		}
	}


}


